package net.yegong.mva.pca;

import java.util.LinkedList;
import java.util.List;

import net.yegong.matrix.DiagonalMatrix;
import net.yegong.matrix.EigenvalueDecomposition;
import net.yegong.matrix.Matrix;
import net.yegong.matrix.RowVector;
import net.yegong.matrix.SymmetricMatrix;
import net.yegong.matrix.Transformer;
import net.yegong.mva.ComputingException;

/**
 * @author cooper
 */
public abstract class AbstractPCACalc implements PCACalc {

	public PCAResult calculate(Matrix m, int factor, boolean normalization, float minCumulative, float errors) throws ComputingException {
		this.errors = errors;
		this.minCumulative = minCumulative;
		this.normalization = normalization;
		// The missing value treatment will release in next generation
		// data = removeMissingValue(data, nRow);
		calc(m, factor);
		PCAResult ret = new PCAResult(eigenvalue, eigenvector, contribution, cumulative, mean, stdev);
		clear();
		return ret;
	}

	private void clear() {
		eigenvalue = null;
		eigenvector = null;
		contribution = null;
		cumulative = null;
		mean = null;
		stdev = null;
	}

	protected abstract void calc(Matrix m, int factor) throws ComputingException;

	protected float minCumulative;
	protected float errors;
	protected DiagonalMatrix eigenvalue;
	protected Matrix eigenvector;
	protected RowVector contribution;
	protected RowVector cumulative;
	protected RowVector mean;
	protected RowVector stdev;
	protected boolean normalization;
	protected float sum_var;

	/**
	 * compute the factor from minimum cumulative contribution
	 * 
	 * @param cov the covariance matrix
	 * @return the factor
	 */
	protected int calcFactor(SymmetricMatrix cov) {
		EigenvalueDecomposition evd = new EigenvalueDecomposition(cov, true);
		evd.setErrorsControl(errors);
		DiagonalMatrix d = evd.getD();
		sortEigen(null, d, d.getRowsCount());
		float tmp = sum_var;
		sum_var = d.trace();
		calcContribution(d);
		sum_var = tmp;
		for (int i = 0; i < cumulative.getColumnsCount(); ++i) {
			if (cumulative.get(i) >= minCumulative - errors)
				return i + 1;
		}
		return cumulative.getColumnsCount();
	}

	protected void setNegativeComponent(Matrix V, DiagonalMatrix D, int factor) {
		if (D != null && V != null) {
			for (int i = 0; i < factor; ++i) {
				float v = D.get(i);
				if (Float.isNaN(v) || Float.isInfinite(v) || v < 0) {
					V.setColumn(i, 0);
				}
			}
		}
	}

	/**
	 * On exit, the matrix A is centralized or normalized
	 * 
	 * @param mat the Matrix.
	 */
	protected void calcMean(Matrix mat) {
		int none_const_columns = mat.getColumnsCount();
		RowVector avg = mat.meanOfColumn();
		RowVector std = mat.stdevOfColumn();
		mean = avg;
		stdev = std;
		for (int i = 0; i < stdev.getColumnsCount(); ++i) {
			if (Float.isNaN(stdev.get(i)) || stdev.get(i) == 0) {
				stdev.set(i, 1);
				--none_const_columns;
			}
		}

		Transformer.minusEachRow(mat, avg);

		if (normalization) {
			Transformer.divideEachRow(mat, std);
			sum_var = none_const_columns;
		}
		else {
			sum_var = std.sumSquare();
		}
	}

	/**
	 * the original matrix is ascend, now we should re-sort it by descend.
	 * 
	 * @param eigenval the eigenvalue to be sorted, if null it will not throw
	 *            NullPointerException
	 * @param eigenvec the eigenvector to be sorted, if null it will not throw
	 *            NullPointerException
	 * @param factor indicate the array dimesion
	 * @param nCol indicate the array dimesion
	 */
	protected void sortEigen(Matrix V, DiagonalMatrix D, int factor) {
		if (D != null) {
			float[] d = D.getArray();
			float temp;
			for (int i = 0, j = factor - 1; i < j; ++i, --j) {
				temp = d[i];
				d[i] = d[j];
				d[j] = temp;
			}
		}
		if (V != null) {
			float[] v = V.getColumnMajorArray();
			int nRow = V.getRowsCount();
			float[] buf = new float[nRow];
			for (int i = 0, j = factor - 1; i < j; ++i, --j) {
				System.arraycopy(v, i * nRow, buf, 0, nRow);
				System.arraycopy(v, j * nRow, v, i * nRow, nRow);
				System.arraycopy(buf, 0, v, j * nRow, nRow);
			}
		}
		return;
	}

	/**
	 * calculate the contribution and the cumulative contribution from
	 * eigenvalue.
	 * 
	 * @param eigenval the eigen value
	 * @param factor the dimesion
	 * @return the list of contribution and cumulative contribution
	 */
	protected void calcContribution(DiagonalMatrix d) {
		float[] eigenval = d.getArray();
		int factor = eigenval.length;
		contribution = new RowVector(factor);
		cumulative = new RowVector(factor);
		float cum = 0.0f;
		for (int i = 0; i < factor; ++i) {
			float val = eigenval[i] / sum_var;
			contribution.set(i, val);
			cumulative.set(i, cum += val);
		}
	}

	/**
	 * @param x remove the missing value rows
	 * @param nRow the number of rows
	 * @return the non-mising value data
	 * @throws ComputingException if every row contains missing value
	 */
	protected float[] removeMissingValue(float[] x, int nRow) throws ComputingException {
		int nCol = x.length / nRow;
		List<Integer> missingRow = new LinkedList<Integer>();
		out: for (int i = 0; i < nRow; ++i) {
			for (int j = 0; j < nCol; ++j) {
				if (Float.isNaN(x[i * nCol + j])) {
					missingRow.add(i);
					continue out;
				}
			}
		}
		int newnRow = nRow - missingRow.size();
		if (newnRow == 0)
			throw new ComputingException("The every row of input table has missing value.");
		if (missingRow.isEmpty()) {
			return x;
		}
		float[] newx = new float[newnRow * nCol];
		for (int i = 0, j = 0; i < nRow; ++i) {
			if (missingRow.get(0) == i) {
				missingRow.remove(0);
				continue;
			}
			System.arraycopy(x, i * nCol, newx, j * nCol, nCol);
			++j;
		}
		return newx;
	}
}
